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1 Introduction 

This chapter presents an overview of the properties of a Bose-Einstein condensate 
(BEC) trapped in a periodic potential. This system has attracted a wide interest 
in the last years, and a few excellent reviews of the field have already appeared in 
the literature (see, for instance, [1, 2, 3] and references therein). For this reason, 
and because of the huge amount of published results, we do not pretend here 
to be comprehensive, but we will be content to provide a flavor of the richness 
of this subject, together with some useful references. On the other hand, there 
are good reasons for our effort. Probably, the most significant is that BEC in 
periodic potentials is a truly interdisciplinary problem, with obvious connections 
with electrons in crystal lattices, polarons and photons in optical fibers. Moreover, 
the BEC experimentalists have reached such a high level of accuracy to create 
in the lab, so to speak, paradigmatic Hamiltonians, which were first introduced 
as idealized theoretical models to study, among others, dynamical instabilities or 
quantum phase transitions. 

The key feature of our problem is that the periodic potential naturally intro- 
duce a spatial discreteness in a nonlinear medium. The periodic potential is gener- 
ally realized with two counterpropagating laser beams [4, 5, 6, 7, 8, 9, 10, 11, 12], 
so as to create an optical lattice (OL). As expected in the mean-field GPE limit, 
the BEC Bogoliubov excitation spectrum has a band structure, in analogy with 
the electronic Bloch bands [13, 14, 15, 16, 17, 18, 19, 20, 21]. When the the power 
of the laser is fairly larger than the chemical potential, the lowest band dynamics 
maps on a discrete nonlinear Schrodinger (DNLS) equation [22]. This was an in- 
teresting remark especially because the DNLS was already widely investigated per 
se by the nonlinear physics community [23, 24, 25] which, indeed, was immediately 
attracted by the new possibilities offered by this system. 

The BEC GPE dynamics in the array can be therefore studied in the frame- 
work of the nonlinear lattice theory [22, 26, 27, 28]. The typical confining potential 
is given by the superposition of an harmonic trap and a periodic potential. For a 
ID OL, the frequency at the bottom of the wells is typically of order of ~ kHz 
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in the OL direction, and the transverse confinement is provided by the magnetic 
potential (characterized by frequencies of order of 100Hz). The axial dynamics of 
a Bose condensate induced by an external potential with cylindrical symmetry in 
the transverse directions can be studied introducing an effective ID CPE equation 
[29]: for BEC in OL, assuming that the Wannier wavefunctions (localized in each 
well) can be expressed in Thomas-Fermi approximation, it has been shown [30] 
that the main effect of the transverse confinement is to modify the degree of non- 
linearity of the DNLS equation, giving raise to a generalized version of the DNLS 
equation. BEC in a periodic potential can allow for the observation of intrinsic 
localized modes (i.e. matter excitations localized on few lattice sites), as well as 
the study of solitons and breathers, possibly also with condensates having a repul- 
sive interatomic interaction. We should also mention that the realization of two- 
and three-dimensional optical lattices [6, 31] opens the possibility to study dis- 
crete/nonlinear effects in higher spatial dimensions. A discussion of the derivation 
of the generalized DNLS equation is presented in Section II. 

In free space, the superflow of a uniform BEC is described by plane waves, 
which becomes energetically unstable in presence of defects when the BEC ve- 
locity is faster than sound, which is the Landau criterion for superfluidity. The 
propagation of sound in a harmonically trapped condensate without OL has been 
observed experimentally [32] and studied theoretically [33, 34, 35, 36]. In the pres- 
ence of a periodic potential, the condensate wavefunction can be expanded in Bloch 
waves, having amplitudes modulated with the periodicity of the OL, which can also 
become energetically unstable when the group velocity is larger than the sound 
velocity [17, 20, 37, 38, 39, 40]. The energetic instability manifests itself with the 
emission of quasi-particles out of a condensate flowing against a small obstacle. 
This happens when the condensate velocity is larger than a critical value, which, in 
the limit of small obstacles is the sound velocity. The interplay between discrete- 
ness and nonlinearity is also crucial for the occurrence of modulational instabilities 
(MI), well known in the theory of nonlinear media. MI are dynamical instabili- 
ties characterized by an exponential growth of arbitrarily small fluctuations of a 
carrier wave, as a result of the interplay between dispersion and nonlinearity. The 
consequences of the modulational instability of the motion of BEC wavepackets in 
OL have been discussed in [41]. A different parametric instability, which will not 
be discussed here, can arise when modulating in time the height of the interwell 
barriers or the strength of the intcrparticle interaction [42, 43]. In Section III we 
will discuss the excitation spectra of a BEC in a periodic potential, while in Sec- 
tion IV we review the occurrence of a discrete modulational instability, comparing 
its effects with those of the Landau instability; a brief discussion of the the prop- 
agation of sound in the OL will be also presented. A discussion of the dynamics 
of BEC wavepackets in OL is given in Section V. 

At last (but not at least) the high laser power available nowadays allows for 
the investigation of low tunneling rates between adjacent wells of the periodic po- 
tential: in these regimes the quantum fluctuations play an important role, and, 
with a strength of the optical potential Vq large enough, it is expected a quantum 
transition from a superfluid phase to a Mott insulator phase. The model used to 
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describe the quantum (beyond Gross-Pitaevskii) properties of ultracold atoms in 
deep optical lattices is the Bose-Hubbard Hamiltonian, which is nothing less that 
the quantized version of the DNLS Hamiltonian. The basic energy scales in the 
Bose-Hubbard model are the tunneling energy K (which decreases by increasing 
Vq) and the charging energy [/ 2 , due to the interaction among particles in the same 
well. Since the role of the quantum fluctuations depends on the ratio U2/K, the 
OL provide an unique way to tune the effective interaction by varying the laser 
power Vq. The phase structure is determined by the two competing terms of the 
Hamiltonian [44] : the interaction energy U2 leads to localization of particles in the 
lattice (Mott phase), while the hopping term K favors superfluidity. The phase 
coherence of different condensates in the supcrfluid phase (and its disappearance 
in the Mott regime) in the array plays a crucial role in the dynamics, and can be 
experimentally studied observing the interference patterns created by the conden- 
sates after turning off the trapping potential. The observation of squeezed number 
states was reported in [45] , while an experimental detection of the Mott-supcrfluid 
transition has been reported for 3D OL in [31]. A systematic study of quantum 
phase transitions in low-dimensional (ID and 2D) OL is presented in [46]. We also 
mention that adding a disordered potential - created e.g. by an optical speckle 
potential [47, 48, 49] or by superimposing laser with different periodicity [50, 51] 
- one expects, in presence of deep optical lattices and for large values of U2/K, a 
Bose glass phase [44]. For space reasons, we will not discuss here the main prop- 
erties of the Bose-Hubbard Hamiltonian and we refer the reader to the chapter 
devoted to beyond Gross-Pitaevskii effects. 

2 Discrete Equations for the Dynamics 

The T = dynamics of a BEC in an external potential V{r) follows the GPE 
[52, 53, 54, 55] 

ih^ t = -^^ + [V + g\ij\ 2 }^, (1) 

where g = AirT^a/m, with a the s-wave scattering length and m the atomic mass. 
The condensate wave function is normalized to the total number of particles N. 
We write the external potential V = Vmt + Vol as the sum the optical lattice 
potential Vol, created by two or more counterpropagating laser beams, and the 
trap potential Vmt, whose form depends on the particular realization of the ex- 
periment. For a ID OL, created by only two counterpropagating laser beams, it is 
VolM = Vl(v, z) cos 2 (2-kx/X), where A = A; aser sin {9/2), A; aser being the wave- 
length of the lasers and 9 the angle between the counterpropagating laser beams. 
The spacing in the lattice is d = A/2 and Vl{v, z) is determined by the transverse 
intensity profile of the (nearly gaussian) laser beams. E.g., in [4], A = 850 nm 
and the 1/e 2 radius of the transverse profile is w 80/ito, an order of magnitude 
larger than the transverse radius of the condensate, so that we can approximate 
the periodic potential by 
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Vol(x) = V cos 2 (k x x) (2) 

where k x — 2tt/X, Vq is the trap depth at the center of the beam, and Vq = sEr 
where Er = h 2 k 2 /2m is the recoil energy. A 2D (3D) OL reads Voh(x,y) — 
V [cos 2 (k x x) + cos 2 (k y y)] (V h(r) = V [cos 2 (k x x) + cos 2 (k y y) + cos 2 (k z z)] ). 
For the ID periodic potential (2) it is useful to write V(r) = Vd(x) + Vl(t), 
where Vd(x) is the x component of the potential Vmt = V x (x) + V y (y) + V z (z). 
Vd has a simple physical meaning: F = — ^f- is the effective force acting on the 
center of mass of a condensate wave packet moving in the periodic potential. 

When the laser power (i.e. Vq) is large enough, we can use a tight-binding 
approximation and decompose the condensate order parameter ip(r,t) as a sum of 
wave functions <&(r — rj) localized in each well of the periodic potential: 

iP(r,t) = J2^(t)$(r-r j ), (3) 

where we denote by j the different wells in the array and ipj (t) = J Nj (t) e l( ^ i W 
is the j-th amplitude. Normalizing to 1 the <£'s, it follows \ ipj \ 2 — N. 
By replacing ansatz (3) in (1), the GPE reduces to a DNLS equation [22]: 

ih lk = ~ K + + + U2 1 ^' |2 (4) 

where the tunneling rate is 

K ^-J dr V ^ ' V ^+i + > ( 5 ) 

the on-site energies are e„ = J dr [|^(V0„) 2 + V"0 2 ] and the nonlinear coefficient 
(which we will suppose equal in each site) is 

U 2 =gNjdr$ 4 n . (6) 

Naturally, if one has a 2D (3D) OL, then the ansatz (3) would lead to a 2D (3D) 
DNLS equation. Equation (4) is the equation of motion tpj — d ®^*^ , where TL is 

the Hamiltonian function 

n = -Kj2 + m+i) + E i i 2 +y 1 ^ |4 ) • (7) 

Both the Hamiltonian TL and the normalization are conserved. 

In the tight-binding ansatz (3) one includes only corrections from the first 
band, which is correct for large Vq. For Vo intermediate is useful to introduce 
contributions from the higher bands, i.e. by considering the ansatz ip(r, t) = 

~ipj,j{i) &f{r ~ r j), where 7 labels the bands: a discussion of the resulting 
discrete vector equation is presented in [28]. We also notice that the DNLS equa- 
tion for just two sites describes the dynamics of BEC in a double well, which 
reduces to the dynamics of a non-rigid pendulum [56, 57, 58]: the dynamical split- 
ting of a BEC into two parts has been experimentally studied in [59, 60], while the 
direct observation of the atomic tunneling in a BEC double well has been reported 
in [61]. 
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2.1 Effects of Transverse Confinement 

The assumption (3) firstly relies on the fact that the interwell barrier Vq is much 
higher than the chemical potentials (e.g., in [9] for Vq ~ bE^ it is fj, <~ O.lVb). A 
second important condition is that the energy of the system should be confined 
within the lowest band. Higher energy bands are not contained in the DNLS equa- 
tion, and become important when the energy is of the order of fru, where u> is the 
harmonic frequency of a single well of the lattice. The effective dimensionality of 
the BEC's trapped in each well can also play a crucial role [62, 30], by modifying 
the degree of nonlinearity of the DNLS equation. In this prospect, the DNLS equa- 
tion can be seen as a zero-order (perturbative) approximation of more complicated 
discrete, nonlinear equations. 

The density profile of each condensate can strongly depend on the number of 
atoms present at a given instant in the same well. This introduce site- and time- 
dependent parameters in the DNLS Eq. (4), modifying, in particular, its effective 
degree of nonlinearity. The tight-binding approximation of nonlinear systems has 
to be generalized as [30] 




with <?j(r; Nj(t)) depending implicitly on time through Nj(i) = \ipj(t)\ 2 . We stress 
here, and discuss again later, that the spatial wavefunctions <£j (which arc con- 
sidered sufficiently localized in each well) can also depend explicitly on time due 
to the excitation of internal modes. For typical experimental setups, however, we 
can consider the adiabatic limit in which the interwell number/phase dynamics is 
much slower that the typical time associated with the excitations of such internal 
modes (and, of course, the cases where such modes are not already present in the 
initial configuration of the system). In this limit, which can be well satisfied in 
experiments, the spatial wavefunctions in Eq. (8) will adiabatically follow the tun- 
neling dynamics and can be approximated with the real wavefunction <Pj(r; Nj(t)). 
A discussion of the validity of the adiabatic approximation is in [30]. 

Replacing the nonlinear tight-binding approximation (8) in the GPE (1) and 
integrating out the spatial degrees of freedom one finds the following discrete 
nonlinear equation (DNL) [30]: 



ih^g = -x + ^*_i) + c.c] ^ + €jii>j + a4 oc i>j 

-[K + X (I 4>j I 2 + I | 2 )]^+i ~[K + X{\ 4>j | 2 + | fPj-i \ 2 M-i. (9) 



In Eq. (9), the "local" chemical potential is the sum of three contributions 



fij depends on the atom number Nj through the condensed wavefunction <Pj. 
The tunneling rates Kjj±i between the adjacent sites j and j ' ± 1 also depend, 
in principle, on the respective populations: expanding the wavefunctions around 





(8) 




(10) 
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an average number of atoms per site, No, and keeping only the zero order term 
<Pj (Nj ) ~ $j(N ) one finds Kjj±i w K, with K given by Eq. (5). The relative 
error committed in this approximation is order of 1CT 4 for typical experimental 
setups. The coefficient \ is given by 

X = -9 J dr S 3 3 S J±1 . (11) 

The on-site energies arising from any external potential superimposed to the OL 
are Cj = J dr Vd £j oc j 2 (ej cx j) when the driving field is harmonic (linear) 
- moreover ej does not depend on the on-site atomic populations. Numerical esti- 
mates show that spatial integrals involving next-nearest-neighbor condensates, as 
well as terms proportional to / dr <P 2 &j±i, can be neglected, but not the terms 
proportional to \. E.g., setting £ = g J dr <P 2 4> 2 ±1 one has - for Vb ~ 2QEr and 
N w 10000 - xN /K ~ lO" 1 and (N a /K ~ 10~ 4 . In a double well potential e.g., 
with height barrier Vb ~ 2n ■ 500Hz and No « 3000, one has x^o/K ~ 1, while 
QNq/K ~ 10~ 3 . For these reasons, one cannot neglect the \ terms in Eq. (9). 
Further studies of a BEC in a double well potential without neglecting terms 
proportional to J dr <P 2 ^|±i are presented in [63, 64]. 

To make Eq. (9) useful, one has to guess the dependence of the localized 
wavefunction 4>j on Nj. It turns out that a reasonable choice is given by sup- 
posing a Thomas- Fermi expression for the <£j's: to be more explicit, let us intro- 
duce the potential V at the bottom of wells, obtained expanding the potential 
V around the minima. At the lowest order V w (m/2)(uj 2 x 2 + &yU 2 + & 2 z 2 )- 
One has to compare the interaction energy with the frequencies w Ii9 , z : we de- 
note by V = 0, 1, 2, 3 the number of spatial dimensions in which one can use the 
Thomas- Fermi approximation. E.g., T> = 3 means that we can approximate <Pj 
with the Thomas-Fermi expression <Pj(r;Nj) oc (jlj — V(r)) where jlj is fixed by 
the normalization condition and depends on Nj - with T> = 2, denoting by (let 
say) y and z the directions in which one can apply the Thomas-Fermi expres- 
sion, one can factorize <Pj = (j)j(x)<p^ F (y, z) with the Thomas-Fermi expression 
(f>T F (y,z;Nj) cx (jlj — (m/2)(u> 2 y 2 + uj 2 z 2 ), jlj yet being determined by the nor- 
malization of the &j. Proceeding along this way, one gets [30] 

n] oc = U a \i,j\ a ; «=2T©- (12) 

The coefficient U a is obtained from Eq.(10) and depends in general on the specific 
trap potential. An estimate for it in a particular setup is given below in Eq.(14). 
The DNLS Eq. (4) is recovered from the DNL Eq. (9) in the case V = (i.e. 
a = 2) and neglecting terms proportional to \. In conclusion the main effect of the 
transverse confinement is to change the degree of nonlinearity and the generalized 
DNLS equation reads 
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where T = [ipj{^* +1 + V'j-i) + c.c] i/jj + {Nj + N j+1 )]i/) j+1 + {Nj + A^_i)]^_!. 

To make the previous result more transparent, let us consider an harmonic trap 
potential Vmt = {m/2)(uj 2 x 2 +u: 2 y 2 +oj 2 z z 2 ) . When the $j does not depend on Nj, 
one has V = and a = 2, as in the standard DNLS equation. However, for deep 
ID lattices, the effective frequencies in the x direction is given by uj x = -\/2Vbfc 2 /m 
and is u) x <~ 10kHz for Vq ~ 5Er, while u) Xi y iZ /2n <~ 100-ffz, for 87 i?6. Then, for a 
number of particles <~ 1000 — 10000 one can use a Thomas- Fermi dependence on 
iVj for the wavefunctions in the y and z directions, but not in the x direction: with 
the previous notation, this means V = 2 and a = 1. This result can be simply 
obtained by factorizing the localized wavefunction <Pj as a product of a gaussian 
<pj having width a (in the x direction) and a Thomas-Fermi cj)j) F (in the y and z 
coordinates): replacing in Eq. (1) and integrating out along the x direction, one 
obtains ej = Qj 2 , where Q = ^muj 2 x {^) 2 , getting the DNL (9) with V = 2 (a = 1) 
and 

Ui = J mLj 2 g/V2irir<j. (14) 



3 Excitation Spectra 

In this Section the Bloch and the Bogoliubov excitation spectra of the system in 
absence of any driving field (Vl = 0) are derived in the tight binding approxima- 
tion. We also present a brief discussion of the comparison with numerical results 
for the excitation spectra of the continuous GPE [20] . 

3.1 Bloch Spectrum 

The Bloch states \P p (r) = e^^&pfr), where *Pp(r) is periodic in the x direction 
with period d, are exact stationary solutions of the Gross-Pitacvskii equation (1). 
The energy per particle e 1 {p) (Bloch energy) and the chemical potential /Lt 7 (p) of 
such solutions form a band structure, so that they can be labeled by the quasi- 
momentum p and the band index 7. 

The generalized DNLS equation (13) describes only the lowest band of the 
spectrum. Exact solutions of the DNL equation are the "plane waves" ipj — 
tpo e l ( fe J-M*)/ft ; where p — Tik/d is the quasi-momentum. Note that the ipj are 
plane waves in the lattice, but do not correspond to plane waves in real space. 
Within the DNL equation framework, the energy per particle e(k) and chemical 
potential /i(k) corresponding to these solutions are found to be [20] 

T) 2 

e(k) = e loc -2(K + 2 X N ) cos (k) = e loc - cos (k) , (15) 

h 2 

H(k) = v loc -2(K + 4 X N ) cos (k) = n l ° c - ~r 9 cos (jfe) , (16) 

where e toc = 2U a N* /2 /(a + 2) and fi loc = ^°% l= ^ Q = d(N a e toc )/dN , with N a = 
\ipo\ 2 the number of atoms per well. In the previous equations we have introduced 



8 



A. Smerzi and A. Trombettoni 



the effective masses m e and m^, to emphasize the low momenta (long wavelength) 
quadratic behaviour of the Bloch energy spectrum and of the chemical potential 
[21]. It turns out that several dynamical properties of the system can be intuitively 
understood in terms of such effective masses. This approach is quite common, for 
instance, in the theory of metals, where m M = m e . However in BEC, because of 
the nonlinearity of the Gross-Pitaevskii equation, the two relevant energies of the 
system, e and p, have the same cos (k) dependence on the quasi-momcntum p, but 
different curvatures. Therefore, m M ^ m e , with 



I if: 2d 2 (K + 2 X N ) 1 _ d 2 p _ 2d 2 (K + 4 x N Q ) 

h 2 



m £ dp 2 



Ti 2 m M dp 2 



It is possible to extend the definition of the effective masses to the full Brillouin 
zone, introducing the quasi- momentum dependent masses m £ (k) = (d 2 p/dp 2 )^ 1 = 
m e / cos(fc) and m l _ l (k) = (d 2 p/dp 2 )^ 1 — m At /cos(fc), where m £ = m £ (0) and 
m M = m M (0). 

Similarly, one can introduce two different group velocities, defined as 

de 1 h dfj, 1 h . 

v e = — = — -sin (k , v fi = — = -sin (k . (18) 

op m e a op m M a 

These two different group velocities are related by [20, 21] — v e + §^N 
with, given Eqs. (17), > v E . The current carried by a Bloch waves with quasi- 
momentum p is p n v e (p), where po is the average particle density; m M , on the other 
hand, plays a crucial role in the Bogoliubov spectrum, which we will discuss below. 

The concept of effective mass, defined as the inverse of the curvature of the 
corresponding spectrum (as that of group velocity, defined as the first derivative) 
can be extended to shallow OL, where the nonlinear tight binding approximation 
breaks down. In this case, the quasi-momcntum dependence of e and p will not be 
simply described by a cosine function, but will still remain periodic in the quasi- 
momentum p. In particular, the value k where m e (p) changes sign (corresponding 
to d 2 e/dp 2 = 0) will be greater than tt/2 and will in general not coincide with the 
momentum where m M changes sign (corresponding to d 2 p,/dp 2 = 0). 

We remark that the Bloch states are not the only stationary solutions of the 
Gross-Pitaevskii equation. Because of nonlinearity, indeed, periodic solitonic so- 
lutions can also appear for a weak enough periodic potential, introducing new 
branches in the excitation spectra [65]. 



3.2 Bogoliubov Spectrum 

In this subsection we study the Bogoliubov spectrum of elementary excitations. 
This describes the energy of small perturbations with quasi-momentum q on top 
of a macroscopically populated state with quasi-momcntum p [stationary solution 
ofEq. (1)]. 

We consider first the case x = 0: in the homogeneous limit (ej — 0), the 
stationary solutions of Eq. (13) are plane waves tpj(t) = ipo exp [i(kj — vt)\, 
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with frequency v given by Tiv — —2Kcos(k) + U\ipo\ a - The stability analysis 
of such states can be carried out by perturbing the carrier wave as ipj(t) = 
(ipo + u{t)e tq ' J + v*(t)e~ tqj ) e »( fe J- I/ *). Retaining only terms proportional to u/tp 
and v/tpo, one gets 

with A = 2K cos(k)-2K cos(fc + q)+{l/2)Ua\^\ a andC = (1/2)[/cm/>o Q/2 ~Vo /2 " 
[20]. From Eq. (19) it follows that the excitation spectrum (i.e., the Bogoliubov 
dispersion relation) for the DNLS with nonlinearity degree a is: 



u±/2K = sin (jfe) sin (g) ± y^cos 2 (Jfe) sin 4 (|) + ^r|Vo| Q cos (k) sin 2 (|) . (20) 

The carrier wave becomes modulationally unstable when the cigcnfrequcncy u> in 
Eq. (20) becomes imaginary: the condition for stability is 

4ifcos 2 (fc)sin 2 (|) +aU\ij)v\ a cos{k) > 0. (21) 

When U is negative (positive), corresponding to negative (positive) scattering 
length, the plane waves with cos (k) < (cos (k) > 0) are stable. When the lhs side 
of Eq. (21) becomes negative, as a consequence of the fact that eigenfrequency w in 
Eq. (20) becomes imaginary, there is an exponential growth of small perturbations 
of the carrier wave: we refer to this instability as the modulational instability. 
In the general case one can repeat the previous stability analysis getting 

[20] 

h 2 sin(fc) sm(q) , „ h 4 cos 2 (fc) sin 4 (o) h 2 N da . 2/ N ,_„. 

R " ±W ±2 V ^^ + ^dk C0S{k)Sm [q) (22) 

with the chemical potential given by u = ji loc — J? m cos (fc) (see Eq. (16)), and 
^loc = u a \ip \ a . For a = 2 (i.e. V = 0) and in the limit x = 0, we recover the 
well known results for the discrete nonlinear Schrodingcr equation [66] . 



4 Landau and Dynamical Instabilities 

From the relation (22), valid for the DNL (13) with x / 0, the small q (large 
wavelength) limit of the Bogoliubov dispersion relation becomes 

hu^ -^-sin(fc) q+\q\J— N cos (k), (23) 
dm^ V m e 8N 

(we assume, for the moment, that ^--§^N cos(fc) > 0). The linear behaviour in 
q indicates that the system supports (low amplitude) sound waves, propagating 
on top of large amplitude traveling waves with velocity 
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. duo 

v s ,± = h-z- 

oq 




(24) 



where the "chemical potential group velocity" has been defined in Eq. (18), and 
the "relative sound velocity" c is defined as 



1 dfi -N Q cos(k). (25) 



The two velocities v Si ± correspond, respectively, to a sound wave propagating in 
the same and in the opposite direction of the large amplitude traveling wave. 

We remark that, contrary to the case of a Galilean invariant system (s = 0), the 
sound velocity depends on the quasi-momentum p. Moreover, v s depends on the 
effective dimensionality of the condensates, since (cf. Eqs. (12) and (16)) ^§^Nq <~ 

a U a N^ 2 . In the limit a = 2, p — > and m s , m M — > m we get the sound velocity 
in the uniform case. 

The system is energetically unstable if there exists an w < 0. In the limit 
s = 0, this corresponds to a group velocity larger than the sound velocity (Landau 
criterion for superfluidity) . When the system has a discrete translational invariance 
(s > 0) the condition for this instability is obtained from the Bogoliubov excitation 
spectrum Eq. (22). Then, we have that the system is not superfluid when u> < 0, 
corresponding to v 2 > c 2 . This result should be compared with the well known 
Landau criterion for an homogeneous system (s = 0), stating that the superfluid 
is energetically unstable when v 2 > c 2 , v = = ^ being the group velocity of 

the condensate, and c = ^i^No the sound velocity. 

There is a further dynamical (modulational) instability mechanism associated 
with the appearance of an imaginary component in the Bogoliubov frequencies, 
which disappears in the absence of interatomic interactions, or in the translational 
invariant limit (if a > 0). The onset of this instability in the tight binding regime, 
coincides with the condition 



c 2 



7T 



<0 cos(fc)<0 => \k\ > - . (26) 



The dynamical instability drives an exponentially fast increase of the amplitude 
of the - initially small - fluctuations of the condensate (while the energetic insta- 
bility should manifest itself in polynomial time [67]). Since the initial phases and 
amplitudes of the fluctuation modes are essentially random, their growth induce a 
strong dephasing of the condensate, and dissipates its translational kinetic energy 
(which is transformed in incoherent collective and single particles excitations) . We 
remark here the different scaling of the energetic and dynamical instability with 
the interatomic interactions. Decreasing the scattering length, the sound velocity 
decreases, and smaller and smaller group velocities can break down the superfluid- 
ity of the system (when a —> 0, the sound velocity c — > 0: in the limit of vanishing 
interactions the condensate is energetically unstable for an arbitrary small group 
velocity). On the other hand, the dynamical modulational instability criterion does 
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not depend on the scattering length. This apparent paradox is simply solved notic- 
ing that the growth time of the unstable modes actually depends on interactions, 
and diverges when the scattering length vanishes (r — > oo when a — * 0). There- 
fore, a noninteracting condensate is always dynamically stable. There is a further 
point to remark: if we consider a condensate moving with an increasing velocity, 
the system always becomes first energetically unstable, then it hits the dynamical 
instability. As a matter of fact, however, in real experiments the energetic instabil- 
ity can grow quite slowly (and at zero temperature only in presence of impurities 
[17]), so that the dominant dephasing mechanism is given by the modulational 
instability. 



5 Wave-Packet Dynamics 

In this Section we review the main properties to the wave-packet dynamics of a 
BEC in an OL, summarizing here the results of a variational approach, previously 
considered in [22, 30]. The approach uses a general variational wavefunction 

il>j = ^)/|l_lJ e 'f&-«+'l^fl ! (27) 

where £(i) and a(t) are, respectively, the center and the width of the wavepacket, 
p(t) and 5(t) their associated momenta and IC(a) a normalization factor (such that 
J2j Nj = N). / is a generic function, even in the variable X = (j — £)/a. For sim- 
plicity, we will confine ourself to an exponential trial wavefunction f(X) = e~ x 
for the standard DNLS equation (4), i.e. a — 2 and x = 0. A discussion of the gen- 
eral case is reported in [41]. The wave packet dynamical evolution can be obtained 
by using the Euler-Lagrange equations for the Lagrangian C = ^ihipjip* — 7i, 
with H. given by Eq. (7). In the following, we rescale the time as t — > ht/2K, 
measuring the energies in units 2K. We also set A = U2/2K and Ej = ej/2K. 
The equations of motion for the variational parameters are [22] 

V = -fjf ; i = sinp • e""j = cosp^ - S^e^ + - 8§£, 

f = 2r8 cosp • e" 7 * (28) 

where r = a 2 , 77 = l/2_T+_T5 2 /8 and the effective potential V is given by V(r, £) — 
K dnE n exp(— 2{n — £) 2 /_T). The pairs £,p and ^, 5 are canonically conjugate 
dynamical variables with respect to the effective Hamiltonian 

H = -cosp • e- ? ' + V(£,r). (29) 

2v ttT 

The effective mass m* = m £ = (since x — 0) is given by ^_ = <LH- — cosp e~ v : 
the quasi-momentum dependence of the effective mass allows a rich variety of 
dynamical regimes. Solitonic solutions with a positive nonlinear parameter A > 0, 
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Fig. 1. Plot of the wave function den- 
sity p n —\ ip„ | 2 at times t = 0,20,40 
with A = 1 in the diffusive regime. Nu- 
merical values: po = 0, 5o = 0, Jo = 
50. The critical value of A is in this 
case A c = 24.8. Solid lines: solutions of 
Eq. (4) with 73 sites; dashed lines: solu- 
tions of variational Eqs. (28). 

for instance, are allowed by a negative effective mass. A regime with a diverging 
effective mass m* — > oo leads to a self-trapping of the wave packet, which has been 
recently experimentally observed [68]. 

In the homogeneous lattice, only the optical potential is present (Vmt = Vd = 
0). Therefore the on-site energies E n , as well as V, are constant. The momentum 
is, of course, conserved and it is equal to the initial value: p(t) = p(0) = po- We will 
discuss here only the case A > 0, in order to make contact with the experiments in 
which 87 Kb atoms with positive scattering length a are used; however, we observe 
that the equations of motion (28) are invariant with respect to the replacement 
A — > —A, po — ► po + it and t — > —t. 

A detailed study of the variational equations of motion is in [69] . Here we quote 
only the main results and we discuss rather the physical implications and the com- 
parison with a full numerical analysis. This comparison is surprisingly successful in 
describing even details of the quite complex dynamical and collisional behaviour. 
Stability phase diagrams for such states are obtained by inspection of the pro- 
file dynamics equations [22]. The parameter A is the ratio between the nonlinear 
coefficient, induced by the interatomic interactions, and the coupling between con- 
densates in neighbour wells: it is the only (geometry dependent) parameter which 
governs the dynamical regimes of the system. When A is small, the wave packet 
spreads out; in the opposite limit, the nonlinearity leads to a localization of the 
wave packet. When cosp < 0, an intermediate regime arises: in this case, the 
effective mass is negative and, for a suitable values of A, a balance can be reached 
between nonlinearity and diffusion. In terms of the variational parameters, this 
means that in the diffusive regime, r — > oo and (if po ^ 0) £ — > oo, with an ef- 
fective mass always finite. On the contrary, in the self-trapped regime, r remains 
finite and the center of mass £ cannot go to oo; furthermore, 1/m* — > 0, meaning 
that r) — > oo and S — > oo. Therefore in this regime there is an energy transfer 
to the internal modes of oscillations, since 6 is the momentum associated to the 
wave packet width: in the full numerical solution of Eq. (4) , this corresponds to a 
breakdown of the wave packet. We note that a nonlinear self-trapping occurs also 
in a two-site model [57, 58, 70]. 

When cospo > 0, the solitonic regime is forbidden and we have only the diffu- 
sive and the self-trapped regimes. In order to show the transition between them, 
let us consider first the case po = 0, in which the center of mass of the wave packet 
does not move (£ = 0). Using as initial values <5o — and Jo, the initial value 
of the Hamiltonian (29) is H — A/2y/nTo — e~ 1//2r °. Since the Hamiltonian is a 
conserved quantity, it is H Q = A/2\f^T- e -i/2r-r5 2 /8_ Therefore —4= - H > 0: 




t=4<) 
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Fig. 2. Plot of the wave function density at times t = 0, 10, 20, 30 with T = 100 in the 
self-trapping region. The numerical values of the remaining parameters are as in Fig. 1. 



when Hq > 0, r have to remain finite and the we have a self-trapped regime in 
which the wave packet remains localized and the nonlinearity forbids the diffusion. 
Vice versa, when Hq <0,r , ^oofori^oo: the wave function spreads out and 
we are in the diffusive regime. The transition occurs at Hq — 0, with 

A c = 2^/^Tae- 1/2r " . (30) 

In Figures 1 and 2 we plot the density | ip n | 2 for different times with A in the 
diffusive region (Fig. 1) and in the self-trapped one (Fig. 2): the solid lines are the 
numerical solutions of Eq. (4) , the dashed lines are the solutions of the variational 
equations (28). As we can see from Fig. 2, the numerical solution of Eq. (4) in the 
self-trapping region loses its gaussian shape [68]. From numerical simulations is 
also seen that the occurrence of the transition between the diffusive and the self- 
trapped regimes does not depend on the chosen initial conditions: what is changing 
is the critical value (30). 

Also when p ^ 0, in which the center of the wave packet moves on the lattice, 
there are two distinct regimes. Hq > 0, i.e., r(t) < r max corresponds to the self- 
trapped regime in which the boson wave packet remains localized around few sites, 
while a diffusive regime occurs when — cosp < Hq < 0. In this case r(t — > oo) — > 
oo and £ w — Hq/ tanpQ = const. The transition between the regimes occurs at 
A c = 2v / 7r7o cospo e~ 1//2/ ". With A > A c , the ratio between the initial value of 
the width cr and the limit width <r max (t — > oo) is given by 

a /a mm = (A-A c )/A. (31) 

We checked the stability of the self-trapping transition also considering different 
initial forms of the wave packet. In Fig. 3 we consider a self-trapped state (yl > A c ): 
the variational prediction is that £ — > and that £ — > const. As time progresses, 
the width increases (and it goes asymptotically to a constant value) and the mo- 
mentum conjugate to the width goes to infinity. The full numerical solution cannot 
go to this state, because the transfer of energy to the internal state breaks down 
the wave packet: when the average position approaches to value predicted from 
the variational analysis (thick line), the wave packet deforms until it breaks. In 
the inset we compare the numerical and the variational average position, where 
this deformation determines a deviation between the two lines. We observe that, 
despite the fact that the variational analysis cannot exactly follow the full dynam- 
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Fig. 3. Density profiles at times t = 
0, 1.25, 2.5, 3.75 (solid lines) and at t = 5 
(dotted) for p = tt/4, r = 100 and A = 
50 (A c = 24.8). The thick line represents 
the asymptotic value predicted from the 
variational analysis. In the inset the varia- 
tional (dashed line) and numerical (solid) 
average position vs. time are plotted. 




Fig. 4. Width (dotted line) and 
average position (solid line) calcu- 
lated numerically for A = A so i and 
po = Stv/4 in a finite array of 73 
sites. 




ics in the self-trapping, it can, however, predict the occurrence of the transition 
and give a fairly accurate estimate of the critical point. 

For cospo < 0, soliton-like structures are present (cf. [23, 24, 25] for more 
references on discrete solitons and intrinsic localized excitations). When cosp < 
the self-trapping condition is given by H >| cospo I and the critical value is 

A c = 2y^Ta | cospo I (1 - e~ 1/2ro ). (32) 

For A < A c , r — ► co, while for A > A c , r remains finite. A soliton solution can be 
determined by imposing P = 5 = 0. One finds [22] 

A sol = 2^JF Q I cospo | e-V 2r °. (33) 

For A = A so i the center of the wave packet moves with a constant velocity £ and 
its width remains essentially constant in time. We observe that for i"b > 1, it is 
A c < A so i . In Fig. 4 we plot the average position and the width for A = A so i . Since 
we are not using periodic boundary conditions, when the wave packet arrives to the 
end of the lattice, it hits a wall and upon rebounding, it regains its original shape. 
For A c < A < A so i, £ — > oo while r(t) oscillates, corresponding to a breather 
solution. When Jo > 1, the breather region extends until Ab re ath > A so i [69]. 

Before concluding, we discuss the relation between the present results and 
the modulational instability. As we discussed, if we consider a small perturbation 
on a plane wave ipj oc e lp ° 3 , stability analysis shows that when cospo < the 
eigenfrcqucncics of the linear modes become imaginary driving an exponential 
growth of small perturbations. When cospo > 0, the plane wave is stable. In the 
present case, we are considering not a plane wave, but a localized wave function; 
therefore, we may expect to find the result previously stated by considering the 
case r » 1. In this case, when cospo > 0, A c — > oo and the self-trapped region 
disappears: this corresponds to the result which there is no modulational instability 
for cospo > 0. To the contrary, when cospo < 0, then A c — > and always the 
system exhibits instability to small perturbation: this means that we have only 
self-trapping, as expected. 

To conclude this Section, we observe that the variational approach can be 
applied also in vertical arrays (i.e., when the gravity is acting) when Bloch oscil- 
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lations occur. Similarly, Bloch oscillations are also possible in horizontal optical 
lattices realized by two countcrpropagating laser beams with a frequency detuning 
varying linearly in time [8, 71]. The DNLS description of the dynamics is confined 
to the first band, and then a complete description of the Bloch oscillations in a 
tilted potential requires the study of the continuous GPE (see [72, 37, 73, 74] and 
references therein). A discussion on the Landau-Zener tunneling is presented in 
the following chapter. In an harmonic trap, one can induce and study dipole os- 
cillations suddenly moving the magnetic potential: if the initial trap displacement 
is smaller than a critical value, it is possible to observe coherent Josephson-like 
oscillations [9]. When the initial displacement is larger than a critical value, the 
modulational instability [41] breaks down the dipole oscillations [75]. 
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